% -----------
% defining the variables and parameters
% -----------
close all;
var y pai r rr z g epsr ym paim rm;

varexo etag etaz etar;

parameters
beta paistar rstar tau kappa phii phipai phiy old_ylead old_ylag old_rpai
new_ylead new_ylag new_rpai
rhoz rhog rhor
lambda elas mu teta corinov gammap gammac;

 
% -----------
% setting the parameter values
% -----------

% Log data density is -349.078233.
%  
% parameters
%            prior mean post. mean conf. interval prior pstdev
% 
%       lambda   0.350   0.4949  0.4082  0.5800 beta  0.1000
%       gammac   0.500   0.5233  0.0961  0.8699 beta  0.2000
%         elas   3.000   2.7855  1.9994  3.5885 gamm  0.5000
%       phipai   0.500   0.3964  0.2038  0.5701 gamm  0.2000
%         phii   0.500   0.8445  0.7873  0.9035 beta  0.2000
%         phiy   0.250   0.4078  0.2297  0.5821 gamm  0.1500
%      paistar   4.000   3.9230  2.3599  5.5709 gamm  2.0000
%        rstar   2.000   1.3806  0.6105  2.1014 gamm  1.0000
%         rhog   0.700   0.6620  0.5540  0.7761 beta  0.1000
%         rhoz   0.700   0.6385  0.4070  0.8829 beta  0.1000
%  
% standard deviation of shocks
%            prior mean post. mean conf. interval prior pstdev
% 
%         etag   0.380   0.3013  0.1777  0.4497 invg  0.2000
%         etaz   1.000   0.9678  0.8235  1.1159 invg  0.5200
%         etar   0.310   0.1774  0.1519  0.1996 invg  0.1600
%  
% correlation of shocks
%            prior mean post. mean conf. interval prior pstdev
% 
%    etaz,etag   0.000   0.5256 -0.1329  0.9802 norm  0.4000


%Steady State Values
lambda=0.4949;
elas=2.7855;
teta=0.75;
rstar=1.3806;
beta=(1+rstar/100)^(-1/4);
paistar=3.9230;
mu=0.2;
gammap=0.5;
gammac=0.5233;
%Policy Parameters
phii=0.8445;
phipai=0.3964;
phiy=0.4078;

%AR(1)
rhog=0.6620;
rhoz=0.6385;

%Standard errors

sdz=0.9678;
sdg=0.3013;
sdr=0.1774;
corinov=0.5256;


% -----------

model(linear);

y=((1-(lambda*elas/((1-lambda)*(1+mu)))*(1+gammac*mu/(1+elas*(1-gammac))))
/(1+gammac-(lambda*elas/((1-lambda)*(1+mu)))*(1+(gammac*(1+2*mu))/(1+elas*(1-gammac)))))*y(+1)
+((gammac*(1-lambda*elas/((1-lambda)*(1+elas*(1-gammac)))))
/(1+gammac-(lambda*elas/((1-lambda)*(1+mu)))*(1+(gammac*(1+2*mu))/(1+elas*(1-gammac)))))*y(-1)
-((1-gammac)/(1+gammac-(lambda*elas/((1-lambda)*(1+mu)))*(1+(gammac*(1+2*mu))/(1+elas*(1-gammac)))))*(r-pai(+1))+g;

pai=(gammap/(1+(1+rstar/100)^(-1/4)*gammap))*pai(-1)
+((1+rstar/100)^(-1/4)/(1+(1+rstar/100)^(-1/4)*gammap))*pai(+1)+
((1-teta)*(1-(1+rstar/100)^(-1/4)*teta)/teta)*((1/(1-gammac))+(elas/(1+mu)))*(y-z)
-((1-teta)*(1-(1+rstar/100)^(-1/4)*teta)/teta)*(gammac/(1-gammac))*(y(-1));


r=phii*r(-1)+(1-phii)*(phipai*pai(+1)+phiy*y)+epsr;

z=rhoz*z(-1)+etaz;

g=rhog*g(-1)+etag;

epsr=etar;

ym=y;

rr=r-pai;

paim=paistar+4*pai;

rm=rstar+paistar+4*r;

end;
  


steady;

check;
shocks;
var etag = sdg^2;
var etaz = sdz^2;
var etar = sdr^2;
var etag,etaz = (corinov*sdz*sdg);
end;

periods=1000000;
stoch_simul(order=1,irf=0,simul);

